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In this review, an overview of structural properties and phase transitions in finite spherical dusty (complex) 
plasma crystals - so-called Yukawa balls - is given. These novel kinds of Wigner crystals can be directly 
analyzed experimentally with video cameras. The experiments clearly reveal a shell structure and allow to 
determine the shell populations, to observe metastable states and transitions between configurations as well as 
phase transitions. The experimental observations of the static properties are well explained by a rather simple 
theoretical model which treats the dust particles as being confined by a parabolic potential and interacting via 
an isotropic Yukawa pair potential. The excitation properties of the Yukawa balls such as normal modes and the 
dynamic behavior, including the time-dependent formation of the crystal requires, in addition, to include the 
effect of friction between the dust particles and the neutral gas. Aside from first-principle molecular dynamics 
and Monte Carlo simulations several analytical approaches are reviewed which include shell models and a 
continuum theory. A summary of recent results and theory-experiment comparisons is given and questions for 
future research activities are outlined. 
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1 Introduction 

Complex plasmas are gas plasmas consisting of electrons, ions, and neutral atoms that additionally contain mi- 
croscopic particles with sizes ranging from lOnm to some 10/im. This state of matter is ubiquitous in space, e.g. 
in the interplanetary medium, in interstellar clouds, in comet tails, and in the ring systems of the giant planets as 
well as in mesospheric noctilucent clouds [1, 2]. At the same time, in microchip manufacturing, avoiding particle 
contamination during the many production steps that involve plasmas is a technological challenge [3,4]. On the 
other hand, the growth, transport, and deposition of nanoparticles is the central goal of many plasma deposition 
techniques, e.g. in the manufacturing of amorphous solar cells [5-8]. 

Since the early 1990s, research in complex plasmas is rapidly evolving. The original name dusty plasmas 
is nowadays often replaced by "complex plasmas" - in analogy to complex fluids, in order to emphasize many 
fvmdamentally different properties of this medium that distinguish it from ordinary gas plasmas. The field of 
complex plasmas is now maturing and the interested reader can find timely monographs, e.g. on dusty plasmas 
in space [9], technological applications [10], or on waves in dusty space plasmas [11]. The microparticles in 
complex plasmas are electrically charged by collection of plasma electrons and ions as well as by photoemission 
or secondary electron emission [1,2]. In laboratory plasmas usually the collection processes dominate and the 
particles attain a high negative charge of a few thousand elementary charges. To some extent, complex plas- 
mas behave like negative ion plasmas. In both cases a heavy, negatively charged species substitutes part of the 
electrons and thus affects the mobility and the contribution to shielding by negative charges. 

Most ordinary plasmas in space and laboratory are weakly coupled, which means that the interaction energy of 
nearest neighbours is much smaller than their thermal energy. This situation is completely different in complex 
plasmas with micrometer sized particles which then become strongly coupled [12]: the negative charge can 
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be so large that the interaction energy of nearest neighbours exceeds their thermal energy by several orders of 
magnitude. A consequence of strong coupling is the formation of liquid or solid phases. In 1934, a liquid to solid 
phase transition had been predicted by Eugene WiGNER for the three-dimensional (3D) electron gas in metals 
at low temperature and low density due to strong Coulomb repulsion [13]. This phase transition became known 
as Wigner crystallization, and the solid phase as Wigner or Coulomb crystals. This kind of crystallization has 
been studied for decades in a varity of systems ranging from electrons trapped on top of liquid helium [14, 15], 
electrons trapped in semiconductor quantum wells [16] and quantum dots [17], strongly coupled radio-frequency 
discharges of dusty plasmas [18], hole crystallization in semiconductors [19] to laser-cooled trapped ions in Paul- 
and Penning traps [20,21]. It is remarkable that these physically completely different systems exhibit very similar 
collective behavior which is due to the governing role of the long-range Coulomb interaction. 

What makes the dust systems particularly attractive is that here collective phenomena, such as crystallization, 
occur at a comparatively large length scale of several millimeters and at room temperature which makes detailed 
experimental investigations rather simple. The formation of ordered solid phases in dusty plasmas was predicted 
by Ikezi [22] and first observed as plasma crystals in the early 1990s [23-25]. Since then crystalline structures 
and dynamical processes have been observed in a variety of confinement geometries ranging from ID to 3D, see 
e.g. Ref. [26] for an overview and additional references. The present review will focus on the 3-dimensional 
spherical dust crystals recently observed by Arp et al. [27]. Their crystalline structure was found to be very 
similar to those observed before in ion crystals in Paul traps and the crystals which have been predicted to 
exist in expanding neutral plasmas [28,29]. Compared to the latter, the advantage of these so called Yukawa 
balls is that the individual particle positions and trajectories can be directly visually observed and traced and 
recorded with video cameras. As a result, very detailed comparisons to theory and computer simulations can 
be performed which go far beyond comparisons in more traditional plasmas. In particular, the comparison is not 
limited to average quantities such as distribution functions on long time scales, but also covers correlation effects, 
fluctuations, and even short-time effects become accessible. 

Thus, the aim of this review is to show, on the example of Yukawa balls, that complex plasmas are an excellent 
object to study in great detail static and dynamic phenomena in finite systems such as ground and metastable 
states, transitions between configurations and phase transitions. These properties have, in the past few years, been 
studied in a very close connection between experiment, theory and computer simulation. For this purpose the 
review is organized as follows: In Section [2] we discuss the experimental realization of the Yukawa balls. Then, 
in Section [3] we give an introduction to the theory and computer simulations and discuss the low-temperature 
structural properties of Yukawa balls. Section [6] is devoted to finite temperature excitations of the crystals, to 
structural transitions and the melting behavior. We conclude this paper with a discussion of the results and an 
outline of further research directions in Sec. |2l 



2 Experimental Realization of Yukawa balls 

To realize and investigate spherical 3D dust clouds two ingredients are indispensable: an isotropic confinement 
to prevent Coulomb explosion of the alike charged particles and powerful 3D diagnostics. For a detailed report 
on the progress in the field of 3D diagnostic tools we refer the reader to [30]. Here, we will limit the discussion 
to dust confinement, because the shape of the confinement potential is a central input parameter for theory and 
simulation. 

To confine dust particles inside a laboratory plasma reactor, the gravitational force has to be balanced. Thus, 
dust confinement is typically achieved in the plasma sheath region above an electrode where strong electric 
fields are present [31]. The trapping of the dust particles in horizontal direction is established by depressions of 
the electrode surface [32] or flat metal rings on the electrode [33]. However, this results in a very anisotropic 
confinement potential. The confinement in vertical direction is much stronger and thus these dust clouds are 
mostly 2D systems which nevertheless can form highly ordered crystals with a hexagonal lattice structure [23-25, 
34]. Further, the supersonic ion flow towards the electrode is focused below each particle [35-39]. In multilayer 
systems, the resulting positive space-charge attracts particles in a lower layer. This process is responsible for 
chain formation observed in all dust clouds which are confined in regions of strong electric fields, i.e. regions 
with strong ion flows. This alignment vanishes if small particles and high gas pressure are used [40^2] . However, 
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Fig. 1 (a) Side view of the discharge arrangement for Yukawa balls. The lower electrode is heated (T < 90°C). The vacuum 
vessel is grounded and kept at room temperature. The dust cloud is confined inside a glass cube where the upper and lower 
side are left open. The inset shows an image of a large dust cloud with 1 cm in diameter, (b) A thin slice at the front side of 
the cloud is illuminated. The particles basically arrange in a hexagonal lattice. After Ref. [27]. 

extended homogeneous 3D plasma crystals cannot be generated this way and investigations of dust dynamics are 
not feasible due to strong damping. 

To produce extended 3D dust clouds, various approaches were followed. Merlino and coworkers confined 
dust in a magnetized anodic plasma [43, 44] and investigated dust acoustic waves. Later, similar experiments 
on dust dynamics were performed in other discharges [45-51]. To produce 3D plasma crystals, a number of 
experiments have been performed under microgravity conditions. These experiments have provided many in- 
teresting observations, e.g. of localized crystalline structures [52], of complex plasma boundaries [53,54], of 
coalescence of complex plasma fluids [55], of transport properties [56,57], and of low-frequency waves and insta- 
bilities [58-61]. However, the most striking observation was the formation of a dust-free zone (void) in the center 
of the discharge [52,62]. It was proposed that the ion drag force is responsible for the formation of voids [63-66], 
and a combination of simulations [67-69], experiments [70-73], and recent ion drag models [74-76] were able 
to verify this. Although it was shown very recently that a void closure can be achieved [77], the formation of 
void-free crystalline dust clouds is still an important issue. Such dust clouds in a well defined confinement would 
allow to extend a large fraction of the successful research on 2D plasma crystals to 3D. Interesting observations 
were reported by Annaratone and coworkers [78, 79] who observed spherical dust clouds with less than 50 
particles in a secondary discharge in front of an adaptive electrode. Unfortunately, these clouds are rather in a 
liquid state and their confinement is not yet understood. 

A different approach to confine void-free dust clouds is motivated by experiments, simulations and theories 
on void formation [62,63,73,80]. Those results showed that a void-free dust confinement at moderate plasma 
densities is hampered by the inherent ion flows towards the discharge boundaries. As a consequence, Arp and 
coworkers [27] have modified the usual setup of an asymmetric capacitively coupled rf -discharge [Fig.|2|a)]. To 
establish a vertical temperature gradient the rf-electrode was heated to T = (60 — 80) °C. The grounded vacuum 
vessel was kept at room temperature. The resulting thermophoretic force should balance gravity for particles 
with a radius < 3 /im, i.e. allow for dust levitation in the bulk plasma. In addition to the experimental setup 
used by Rothermel [73], Arp et al. proposed to place a glass box with square cross section on the electrode 
that was left open at the top and at the bottom. For low discharge power (P < 10 W), they discovered that dust, 
which was dropped into the glass box, formed spherical void-free clouds [27]. A picture of a laser illuminated 
dust cloud (Yukawa ball) consisting of roughly 10, 000 dust particles is shown in the inset of Fig. |2|a). Already 
if only the front of such a dust cloud was illuminated with a laser sheet, interesting structural properties were 
observed [Fig.|2|b)]. The particle arrangement was static and particles on the surface of Yukawa balls seemed to 
have preferably 5 or 6 nearest neighbors. Moreover, the particles did not form vertical chains, i.e. first evidence 
was found that the Yukawa balls were in a solid or even crystalline state. 

To achieve a quantitative description of the dust confinement, fluid simulations with the SIGLO-2D code 
[81,82] in combination with additional experiments were performed [83]. The results are summarized in Fig.|2] 
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Fig. 2 Experiments on the confinement of Yukawa balls: (a) Superposition of multiple video frames. As soon as the discharge 
is switched off, the trapped particles fall down. Their motion is only affected by gravity and the thermophoretic force, (b) 
Vertical component of the experimentally determined thermophoretic force field. In 'red' regions the thermophoretic force 
exceeds gravity, in 'blue' regions gravity is dominant. The solid line shows where both forces balance, (c) Trap potential 
obtained from PIV measurements and fluid simulations, (d) A horizontal section through the trap center reveals an almost 
parabolic confinement. FromRef. [83]. 



It was found that thermophoresis alone was not sufficient to explain the confinement [Fig. [2]; a)]. By means of 
particle imaging velocimetry (PIV) experiments the thermophoretic force field was measured [Fig.|2|b)]. It was 
found that thermophoresis contributed about 70 percent to the vertical confinement and that the radial confinement 
and the remaining 30 percent in vertical direction were provided by plasma induced forces, namely electrostatic 
fields due to surface charges on the glass [83]. Further, the fluid simulations showed that the Yukawa balls were 
confined in a region where the plasma production was negligible. Due to the low plasma density (n w 10^'^m^'^) 
and the weak electric fields the ion drifts were significantly below the ion sound speed. The resulting ion drag 
force was about two orders of magnitude smaller than the electrostatic confinement force and could, therefore, 
be neglected and chain formation is not expected. The confinement provided by the combination of gravity, 
thermophoresis and electric fields is plotted in Fig. [SJc) and is in very good agreement with the experimental 
observations [Fig.|2|a)]. It yields the correct levitation height, its depth is sufficient to trap particles with a few 
thousand elementary charges and even its asymmetry was observable for huge clouds in the experiment [83]. 

For further structural analysis and especially for comparison with simulations as well as with other strongly 
coupled systems, it is important to note that the trapping potential is essentially isotropic and parabolic for dust 
clouds with a radius of less than 2 mm, i.e. N < 1000 particles. Furthermore, for the trapping process an 
interaction of dust and plasma, as proposed by TOTSUJI [84, 85], is not required. The particle trap is provided by 
external forces only and thus this trap geometry is closely related to those of Penning and Paul traps [86]. This 
allows for a direct comparison of the structural properties of Yukawa balls and e.g. trapped laser-cooled ions. 
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3 Theoretical Description of Yukawa balls 

A theoretical description of complex plasmas is very complicated, due to the large number of different particle 
species, which include electrons, ions, neutrals and dust particles, and their huge difference in mass. Thus the 
length and time scales of the motion of the different particles differ by many orders of magnitude, making a 
selfconsistent modeling unfeasible. The most advanced computer models use particle in cell (PIC) simulations 
with a full account of the time -dependence of the discharge. [87]. However, this currently allows to treat only a 
very small number of dust particles in the plasma environment. A simplified model which treats the streaming 
electrons and ions within linear response and computes the dynamically screened potential of the dust particles 
was developed by Joyce and Lampe, see [88] and references therein. Both simulations yield an anisotropic 
and non-monotonic (wake) potential around a dust grain if the ion streaming velocity exceeds the sound speed 
which is normally the case in the plasma sheath. At the same time, simulations for the conditions of Yukawa 
balls which are produced in the plasma bulk revealed that there the streaming velocities are much lower and the 
dust-dust interaction is nearly isotropic and statically screened. It is, therefore, well approximated by a Yukawa 
potential. These findings justify the use of a much simpler theoretical model which is discussed in the following 
sections. 

3.1 Theoretical Model 
A standard model to describe the Yukawa balls and their dynamics is based on the Hamiltonian 

JV 2 ^ ^ ^2 

i7(ri...,p„...r^,...p^,i) = ^|^+5:|r,^ + 5:C/-*(r.,t) + ^^.e— ^ (1) 

i—l i—1 i—1 i^j ^-^ 

containing kinetic energy, a confinement and the interaction energy. The isotropic and parabolic confinement is 
based on the experimental results, cf. Sec. |2] q is the particle charge and a the confinement strength. For both 
values the experiment sets rather narrow boundaries [83]. Finally, e is the background dielectric constant and 
is the distance between particles i and j. The effect of the surrounding plasma on the dust is condensed in the 
static screening parameter k which will be treated as a parameter below. A second important effect of the medium 
on the dust particles can be described by a friction force. It is beyond the hamiltonian model but can be included 
in the equations of motion, see Sec. 3.2 In writing the formula ([T]i we assumed that all particles have identical 



mass and charge which is a good approximation of the experimental conditions. In fact, commercially available 
dust particles have a diameter which fluctuates by less than one percent, and the remaining charge variations were 
shown to have a negligable effect on the structure of Yukawa balls [89]. Thus, the present model is expected to 
be well appropriate to reproduce the experimental conditions. 

In Eq. (0 we have also included a possibly time-dependent potential V^^^ which describes external manipula- 
tions of the dust particles. The Hamiltonian ([T} is brought to a simpler dimensionless form suitable for numerical 
analysis by introducing units of length and energy, tq = {q^ /IGireaY^^ and Eq = {aq* /32TT'^e^y^^, where vq 
corresponds to the ground state distance of two particles in the trap in the absence of screening and Eq is the 
corresponding pair interaction energy. As the third scale parameter we introduce the unit of time Tq which is the 
inverse of the trap frequency Tq"^ = lu = \J ajm. 



3.2 Computer simulations 

The model ([T]i is well suited for accurate computer simulations. Since the dust particles consitute a classical 
system, exact simulation methods are available which include molecular dynamics (MD), Langevin dynamics 
(LMD) and Metropolis Monte Carlo (MC) which have been successfully applied in recent years. Since most of 
the theoretical results in this paper are based on these methods we give a brief discussion. 
Molecular dynamics solves Newton's equations corresponding to the hamiltonian ([T]), 

p,(i) = -ar,- Vf/'='^*(r„i)+^-^(l + Kr,,)r,, •e"^'--, i = l,...N, (2) 

r,(0) - r°, pM=pI 
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where the last line contains the initial conditions. The coupled system (j2|i of ordinary differential equations is 
efficiently solved for up to several thousand particles by standard techniques such as Runge-Kutta methods. The 
present version gives a solution in the microcanonical ensemble where the total energy H is conserved. The 
initial conditions are chosen as close as possible to the experimental situation. This scheme has been used to 
obtain the ground and metastable states of the Yukawa balls, see Sec. |4.1| 

To reproduce the equilibrium properties of Yukawa balls which are characterized by a given temperature it is 
more appropriate to perform simulations in the canonical ensemble. This can be realized by Langevin MD where, 
in addition to the total force on particle i [r.h.s. of Eq. Q] a friction force describing the momentum loss of 
the particles to the neutral gas and a random fluctuating force modeling the collisions of the particles with the 
plasma medium (Brownian motion), 

p,(i) = ¥,{t)~vmr,+y,{t) i^l,...N, (3) 
(y*) = 0, {y^a{t)yji3{t')) =2D5^.j5ajj5{t-t'), a,(3 = x,y,z, 

to be supplemented by the same initial conditions as in Eq. Choosing the noise intensity according to the 
Einstein relation, D = mkBTv, the system relaxes towards a Maxwellian velocity distribution with temperature 
T - tpyically the temperature of the neutral gas. This method has been used to obtain the probability of metastable 



states, cf. Sec. 6.2 and for nonequilibrium dynamical simulations of Yukawa balls [90]. 

Finally, a first principle method to study the thermodynamic equilibrium properties of the dust particles is 
Monte Carlo. In the canonical ensemble where the dust temperature Td and the particle number are assumed to 
be fixed, the probability of a certain configuration s of the Yukawa ball with the coordinates R'* = {r^, . . . r^} 
is defined by 

P(R^; Td, TV, a) - ie-«o(R-^)/'c«T,^ j ^3^^ _ d^NPi^'; T^, TV; a) = 1 (4) 

where Z is the normalization constant (partition function). Note that in contrast to thermodynamics of macro- 
scopic systems, here the volume is not fixed but its role is taken over by the confinement strength. The exponent 
iJo of the probability distribution P is the hamiltonian ([l]) with the momenta put equal to zero [the momentum 
and coordinate distributions factorize, and the former is given by a Maxwellian and does not influence the par- 
ticle configurations]. Also, the effect of dissipation is not included in this description. The equilibrium state 
of a Yukawa ball, is obtained by maximizing the probability P which is efficiently done using the Metropolis 
algorithm, e.g. [91]. MC methods have been applied to find the ground and metastable states, see Sec. |4.1[|6.2| 
and finite temperature effects and melting points, Sec.|6.4| 



4 Structure of Yukawa balls 

Let us start the analysis of Yukawa balls by considering their ground state properties. They are obtained in 
the zero-temperature limit where all particles settle in the minima of the total potential energy Hq given by 
the sum of confinement and pair interactions. While zero temperature seems to contradict the experimental 
situation - the experiments are performed at room temperature - this is in fact a very good approximation. The 
reason is the very high charge collected on the grain surface leading to a potential energy per particle which 
exceeds the kinetic energy by a factor of 500 or more. Therefore, kinetic energy has a negligible effect on the 
stationary properties. In most cases the potential energy has several minima. In this case, the deepest minimum 
corresponds to the classical ground state whereas the others correspond to metastable states. Due to the spherical 
symmetry of the hamiltonian i/o the low temperature stationary states have the same symmetry: they consist 
of spherically symmetric concentric shells with a very small width which is in excellent agreement with the 
experiment. Different states differ by the radii and populations Nk of the shells. The particle arrangement within 
the shells is highly symmetric with a combination of five-fold and six-fold symmetry (particles dominantly have 
five or six nearest neighbors) [27,92]. 



4.1 Ground and metastable states 

Much work has been devoted to first-principle simulations of spherical Coulomb and Yukawa crystals. This is 
efficiently done by MD simulations by solving Eqs. Q, combined with "simulated annealing". There one starts 
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with a random initial configuration and reduces the velocity of all particles proportionally in small steps, until 
they have zero kinetic energy, [92,93]. Alternatively, one can use LMD simulations, Eq. (j3]l, with a very small 
damping coefficient v, e.g. [90] or MC. The ground state for the case of pure Coulomb interaction was reported 
many years ago by Hasse et al. [94] where also earlier references are given, for larger crystals, see also [95]. 
These data correctly describe the general trend previously observed in spherical ion crystals. However, with 
the first experimental results on Yukawa balls in dusty plasmas very detailed comparisons of the precise shell 
structure became possible. Then it turned out that most data of Ref. [94] refered to metastable states rather than 
to the ground state. The first accurate data for spherical Coulomb crystals were obtained by LUDWIG et al. [92] 
who performed MD with simulated annealing. Extensive tables for iV < 160 and some larger clusters were 
presented in Ref. [96] and the results were later confirmed by other groups, e.g. [97,98]. 
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Fig. 3 (color online) Stationary states observed in the LMD simulations for A'^ = 31, k = 1.4 for cooUng to a minimum 
kinetic energy {E^in') = 10^**. The runs are sorted by the energy or the final state, see also tablejl] For slow cooling (black 
bars, V = 0.05) one can clearly distinguish different states, cf. the horizontal lines. The length of the bold lines is proportional 
to the occurrence frequency of a state in a total of 5000 runs. In the case of strong friction (red, dashed line, v = 5.3) the 
particles often lose their kinetic energy before they can settle into the equilibrium positions and the fine structure cannot be 
resolved. From Ref. [90] 



AE/N 


config. 


AE/N 


config. 


3.030266 


(27,4) 
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0.000009 


(27,4) 
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(26,5) 


0.000656 


(25,6) 


0.000372 


(26,5) 


0.000669 


(25,6) 



Table 1 Energy difference between metastable states and the ground state (the ground state and its energy is given by bold 
numbers) as seen in Fig. [3] States with the same shell configuration but different energy differ only by the arrangement of the 
particles on the same shell (fine structure) [92]. 



First structural data for Yukawa balls have only recently been obtained, e.g. [99, 100]. The first systematic 
analysis was performed in Ref. [101] for small crystals with < 60 for a range of screening parameters of 
< K < 5.0. An example of the energetically lowest stationary states is given in Fig. |3]and in table [T] for a 
cluster with = 31 and k = 1.4. The table also contains the shell configuration of the states. As one can 
see in some cases there exist several states with the same configuration which differ only in the arrangement 
of the particles on the shells. This "fine structure" was first observed in Ref. [92] and revealed by a Voronoi 
analysis [102]. It is remarkable how small the energetic spacing of these states is which underlines the high 
numerical accuracy needed to find the ground states. Since no method is able to reveal the ground state with one 
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hundred percent reliability it is necessary to perform a large series of repeated calculations, as shown in Fig. |3] 
where 5000 runs were performed. Further, the results sensitively depend on the cooling rate in the annealing or, 
equivalently, on the magnitude of the friction in an LMD simulation. The ground state and the metastable states 
can only be resolved for sufficiently small damping, i.e. for slow cooling [90]. Further note that the number of 
metastable states increases exponentially with N which makes it very difficult to find the ground state for clusters 
with more than about 1000 particles. 




Fig. 4 (color online) Experimental (symbols) and theoretical (lines) shell population Ns of the ground states versus system 
size. The MD results are obtained for different screening parameters k and show that the particles redistribute towards inner 
shells with increased k. From Ref. [99] 



4.2 Effect of screening 

A comparison of experimental and theoretical cluster configurations is shown in Fig.|4]for particle numbers in the 
range of < 500. In addition to the measured configurations of 43 clusters MD simulations of the ground states 
of the same particle numbers and with different values of k were performed [99], using simulated annealing, as 
described in Sec. |4.1| above. Comparison of the results shows that the model of isotropic static screening in fact 
reproduces the measured shell occupations very well. Furthermore, it is obvious that the measurements caimot 
be explained by pure Coulomb interaction. Since k is the only free parameter in the simulations within the model 
([T]l, comparison to the measurements allows for reliable predictions of the screening strengths in the experiments 
which is not directly accessible. The result of k « 0.6, corresponding to A/jf — 1.5, gives a value which agrees 
well with other estimates and the accuracy was estimated to be within 50 percent [99]. Sources of errors such as 
charge fluctuations and finite temperature were found to be small [89] and also the occurence of metastable state 



in the experiments, see Sec. 6.2 has a negligible influence on this result. 

For completeness, we note that Ref. [99] also provided detailed comparisons of the experimental and theoret- 
ical shell radii Rs- It turned out that the radii are almost equidistant and increase approximately as N^^^, as one 
would expect in bulk matter. Most importantly, while the absolute value of the radii decreases with screening 
(due to the reduced interparticle repulsion), no effect of screening is observed for the ratio of Rg to the mean 
interparticle distance within a shell, for fixed N. Since at the same time, the particle number on the outer shells 
decreases with k, this has important consequences for the mean density distribution within the Yukawa balls 
which will be considered more in detail in Sec. 15.21 

There are two general effects of screening on the structure of Yukawa balls: first, increase of screening leads 
to cluster compression and, second, to increased occupation of inner shells. However, this is only the average 
trend, and noteble deviations have been found. A detailed numerical analysis of small Yukawa balls in a broad 
range of screening parameters has revealed a number of interesting anomalies [101], the shell populations are 
shown in Fig. |5] In general the particle number on the inner shell increases one by one, when k is increased. 
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Fig. 5 (color online) Ground states of small Yukawa balls vs. screening parameter. The white numbers denote the number 
of particles on the inner shell(s). The black circles indicate anomalies of the 1st kind. The white circles indicate the end of the 
screening range, where anomalies of the 2nd kind appear. The ground states for a screening parameter k — 20.0 are plotted 
above the diagram. The cyan bar for N — Aiat k. — 20.0 refers to a ground state of (11, 1) in the center region; it is the only 
time this configuration is part of a ground state. The small squares just below k — 20.0 indicate anomalies of the 3rd kind, 
where a ground state configuration reappears with increased screening. From Ref. [101] 



However, in several cases, two particles move inward simultaneously, see for example N — 30 at k sa 1.5 
where the configuration changes according to (26,4) (24,6). These correlated transitions have been called 
"anomaly of the 1st kind" and are observed 18 times. A second anomaly is observed when N is increased by 
one at constant k. The normal trend is an increase of the inner or outer shell population by one. However, in 
some cases the inner shell population decreases by one whereas the outer shell population increases by two. This 
has been observed for Coulomb balls in a single case = 59 ^ 60 [103] when the configuration changes 
according to (46, 12, 1) (48, 12). In the case of Yukawa balls this "2nd anomaly" is observed frequently, see 
e.g. = 11 ^ 12 for 2 < k < 4 where the configuration changes according to (10, 1) (12,0). Finally, 
there are five cases (A^ = 35, 36, 37, 39, 54) of reentrant shell configuration changes where for a fixed A^ and 
increasing k a particle moves from an inner to an outer shell restoring a configuration which existed before, at a 
smaller value of k. The reason for these non-trivial configuration changes are special symmetries of the particle 
arrangement which allow to reduce the total energy of the Yukawa ball. A table containing a complete list of the 
anomaUes is given in Ref. [101]. 



5 Analytical Models for the structure of Yukawa balls 

While computer simulations provide virtually exact structural data for a certain hamiltonian a lot of additional 
information and physical insight can be gained from analytical models. Appropriate simplification of the system 
allows to observe the dominant trends and effects. Here we summarize the results of two such models - shell 
models and statistical models based on a fluid-type continuum description. We start by considering the Coulomb 
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limit which is simpler - the ground state configurations of small clusters with iV < 8 can be found analytically 
[104] - and then proceed to Yukawa interaction. 



5.1 Shell models 



The first idea of a simplified picture of the stationary particle configuration is to use a shell model, where the 
particles are confined to concentric shells with zero width. The first shell model was derived by AviLOV and 
Hasse [94] and was improved by TSURUTA and ICHIMARU [103] who included correlation effects. Kraeft 
and Bonitz further improved this model and presented detailed comparisons with MD simulations [105]. They 
gave the following form of the total energy per particle for a system of L homogeneously charged concentric 
spherical shells 



N{ZeY/ro ~ ^ TV^J^ 

/J,— 1 



1 



10 



.iV2/3 



(5) 



where (^ = 0,1 accounting for a particle in the trap center and is the radius of the shell /i in units of tq. 
Here, the first term (proportional to N^) is the surface energy of a spherical capacitor of radius a, containing N^, 
charges, Esurf{N) = {N{N — l)e^)/(2a). The contribution from the sum over fi accounts for the electrostatic 
interaction of the shell ly with all inner shells. The xl contribution describes the confinement energy and the last 
term takes into account a (possible) compensating background and is not relevant for the structure of the Coulomb 
ball. This model, without the term proportional to e, can be rigorously derived from a mean field theory [106] 



which is discussed in Sec. 5.2 and was given by Hasse et al. [94] although they missed the —1 correction to the 
surface energy Esurf- 

The cluster configuration can now be derived simply by an optimization procedure searching for those shell 
populations and radii which minimize the total energy (|5]l which is much simpler than to solve the exact problem. 
While this yields the correct qualitative trend, however, with e = there are large quantitative deviations. Tsu- 
RUTA and ICHIMARU corrected the Hasse model by introducing a correlation correction (they used e = 1) which 
reflects the discreteness of the particles and takes into account that the area around a particle an a shell cannot be 
occupied by others. Using e = 1 allows to reproduce the shell populations within 5 percent [105]. However, the 
exact expression for the correlation corrections is unknown, and there is no obvious reason why this coefficient 
should be equal to one. In fact, a further improvement of the quality of the model is achieved by leaving e as an 
additional free parameter that is optimized to minimze the energy. This has been done in Ref [105] and allows to 
reduce the deviations from the exact ground states to (1 . . . 2)%. The resulting values for the correlation param- 
eters were slightly above one and converged to e = 1.104 for large N , a result which could be recently derived 
from the Thomson model by ClOSLOWSKl [108, 109]. 

For Yukawa balls the situation is more complex. A mean field-type model (e = 0) was recently obtained by 
TOTSUJI et al. [84] and derived from continuum theory by Henning et al. [106], 



E^odei {N; k) = J2n^{^RI + Q 



V KRt. 2 



-kR^ 



Ri/ 



(6) 



where, compared to Eq. (j5]l, the background term (last term) has been dropped and we multiplied with the energy 
unit (denominator on the l.h.s.). One readily confirms by performing the limit k — > that this result includes 
the Coulomb case. The model (|6]l differs from Refs [84, 106] by the additional correlation corrections which 
generalize the Coulomb expression (j5]l. This model was used in Ref. [107] and optimized to minimize the 
total energy. The detailed comparison to exact simulation data showed that, for Yukawa balls, the mean field 
model (neglecting the e terms) performs rather poorly, but satisfactory agreement with an accuracy of the shell 
populations within 5% can be achieved only by including correlation corrections and using different e-values for 
different shells. Yet a systematic derivation of the correlation corrections remains open. 
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5.2 Continuum theory approach. Radial density profile of Yukawa balls 

An entirely different approach to the many-particle behavior of the Yukawa balls has been developed by Kraeft 
and BONITZ [110] and by Henning et al. in Refs. [106, 111]. It is based on a fluid-like statistical theory where 
the particles are treated as a continuum. This may seem far from reality but one can hope to describe several 
aspects of the Yukawa ball structure such as spatially averaged properties and achieve deeper physical insight. 
One such question is the mean radial density profile of the Yukawa balls. It is well known that, in a spherically 
symmetric parabolic potential, particles interacting via the Coulomb potential establish a radially constant-density 
profile (this follows from the fact that the potential profile inside a homgeneously charged sphere is a parabolic 
function of the distance from the center). However, in the case of Yukawa balls the pair interaction is screened, 
which leads to enhanced populations of inner shells with increaed k while the shell radii (in units of the mean 



interparticle distance within the shells) does not change with k, cf. Sec. 4.2 From this we may expect that radial 
distribution of matter in the Yukawa balls will be inhomogeneous and the density has to decay towards the edge. 
This has indeed been observed experimentally [1 12, 1 13] although the trend is weak, see Fig. [9] 

To analzye this question quantitatively we consider the model of a finite one-component plasma containing N 
identical particles described by the Hamiltonian ([TJ. The classical ground state energy follows from neglecting 
in Eq. ([TJ all particle momenta and can be written as a functional of the density [106] 

(7) 



Ho[n] = E[ji] = J d^rM(r), 
with the potential energy density 

N~ 1 



u{r) = n{r) < <i>(r) + 



2N 



(fr2 n(r2) 



47re|r — r2 



-Klr-ral . , 



(8) 



where the first term on the right is the mean-field contribution with the confinement potential <&, and Ucorr 
denotes the density of the correlation energy. The ground state density profile can be obtained from minimizing 

the total energy E, i.e. from the solution of the variational problem = 5E[n]/ 5n{r) under the constraint 
J (Prn{r) = N =const.This problem can be solved for a general anisotropic confinement potential ^(r) in a 



mean-field approximation, i.e., Ucorr = 0, with the result [106]: 4:TTq'^n{r){N — l)/N 
where /i is a Lagrange multiplier (chemical potential) assuring the normalization. 
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Fig. 6 Radial density profiles of a 3-D plasma of A'^ = 1000 particles calculated with the mean-field model (solid lines) 
and with the LDA including correlation contributions (dashed lines) for Coulomb and Yukawa interaction with four different 
screening parameters: kcLc — 1, 2, 3, and 5. Averaged shell densities of molecular dynamics simulations of a plasma crystal 
for the same parameters are shown by the symbols, dc is the two-particle equilibrium distance ro, cf. Sec. [3] 
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and the density drops to zero at a finite radius R{N, k) which follows from the normalization 

-15 — (iV - 1) - 15 — k(7V - l)R + 15R^ + 15kR^ + Gk^R^ + k^R^ = 0. (10) 
a a 

This equation has four complex and two real solutions, only one of which is non-negative, and thus constitutes the 
unique proper result entering equation (|9]l. This density is equivalent to a local force balance [106, 114]. Figure 
|6] shows the density profiles derived from the mean field (MF) result (j9|l (solid lines) and compares them to the 
averaged shell densities of MD simulations (symbols). In the case of weak screening (a) the MF result agrees 
perfectly with the exact MD results, while for increased screening discrepancies show up, especially in the trap 
center. 

The origin of these deviations is the neglect of correlations. An extension of the previous analysis to include 
correlation effects is difficult as the form of the pair correlation function is unknown. One way around this 
problem is to apply a local density approximation (LDA) [111, 115]. There one replaces the nonlocal terms 
within the energy density at point r by local expressions using the known energy density of the homogeneous 
system [85] 



-0.375Kn//^ + 7.4* IQ- 



(11) 



where no is the corresponding density of the homogeneous system. Using this expression in equation (j8]l, leads 
to an improved equation for the density profile which has to be solved numerically. The resulting density profiles 
for a cluster with 1000 particles are included in Fig. |6] In the limit of strong screening, where the MF result 
failed to describe the simulation results, the LDA approximation with correlations is very accurate, see Fig.|6j). 
In contrast, for weak and moderate screening, kj-q > 2, the LDA behaves very poorly. This is due to the long 
range of the interaction which cannot be accounted for within a local approximation. Thus, the mean-field model 
together with the presented local density approximation complement one another in the description of strongly 
correlated spatially confined one-component plasmas providing a complete theoretical description of the average 
ground state properties of Yukawa balls. An extension of these results to finite temperatures has been outlined in 
Ref. [116] 



6 Structural transitions 

Up to this point we concentrated on the ground state properties of Yukawa balls, but as was shown, in Fig.|3]and 
table[T]also other stationary states - metastable states are frequently resolved in the simulations. Their energy is 
higher than that of the ground state but often only marginally. It may, therefore, be expected that these metastable 
states might be observable in the experiment and also transitions between stationary states might occur. We first 
analyze these questions theoretically and then compare with recent experimental results. 



6.1 Coexistence of stationary states. Potential barriers 

As we have seen in Fig.|5] upon increase of k the ground state configurations frequently change. This requires 
that, around some critical value Kcr, both, the old and new configuration are stationary states, corresponding to 
minima of the total potential energy Hq. When k departs from Kcr one of the configurations (the metastable state) 
will eventually disappear. This is illustrated for iV = 31 in the left part of Fig.|7] The right (left) arrows indicate 
the range of k values where various metastable states appear (vanish). 

If stationary states co-exist one could expect to observe transitions between them which may occur sponta- 
neously, driven by thermal fluctuations. However, this does not only depend on the local minima but also on 
the potential landscape between them, in particular on the potential barrier. The radial potential barriers can be 
"measured" in a computer simulation by moving, e.g. an outer particle to the inner shell, while allowing all other 
particles to relax. The results of such a MC simulation are shown in the right part of Fig.|7] They indicate that, 
at the chosen screening, transitions (6,25) (5,26) and (5,26) (4,27) occur practically spontaneously 
whereas the opposite transitions are hampered by a significantly higher barrier. This asymmetry his even more 
pronounced for the transition between the configurations (4,27) and (3,28). This is nicely confirmed by the 
experiments [113] where a slow transition (4, 27) (5, 26) could be observed followed by a rapid return to the 
ground state (4, 27), see Fig.|9] 



Copyright line will be provided by the publisher 



cpp header will be provided by the publisher 



13 




Fig. 7 (colour online) Left: Energy per particle of a Yukawa ball with 31 particles vs. screenings strength. The ground 
state configuration is (27,4) for k < 1.56, which is replaced by (26,5) in the inverval 1.56 < k < 1.61, followed by 
(25, 6) for larger k. Inset shows the energy differences of the configurations. The vertical dashed lines with the arrows denote 
the beginning [blue for (26,5) and green for (25,6)] and the end [red for (27,4)] of occurance of these configurations in the 
simulations. Right: Potential energy landscape measured by a MC simulation where one particle is moved, while all other 
particles are allowed to relaxate into a new minimum energy configuration (A'^ = 31, k = 0.8). Curves show the potential 
barriers for three transitions (27, 4) (26, 5) (sohd line), (26, 5) (25, 6) (dashed) and (27, 4) (28, 3). The peak in 
the solid curve around 0.75 < ro < 0.9 occurs when the number of MC steps is reduced simulating incomplete relaxation. 
Interestingly this peak matches the return point in the experiment, cf. Fig.|9](a). 

6.2 Probability of occurance of metastable states 

In the case of coexisting stationary states the question arises with what probability each of them will be observed 
when the Yukawa balls are formed. This subtle question could in fact be answered experimentally by performing 
a sequence of repeated measurements where the confinement was rapidly turned off and on again such that 
the plasma parameters were not altered [112, 113]. The statistics was sufficient to extract reUable occurrence 
frequencies of different states which are shown in the left part of Fig. [8]by the horizontal solid straight lines (the 
dashed lines indicate the statistical eiTor). The measurements proved the existence of three stationary states. Most 
suprisingly, the metastable configuration (26, 5) was found to occure with an almost two times higher probability 
than the ground state (27,4) [(62 ± 13)% compared to (35 ± 10)%]. We will now analyze these observations 
theoretically. 

The occurence probability of stationary states can be described analytically by a harmonic approximation of 
the partition function [90]. There, the potential energy of a given state s (ground or metastable state) is expanded 
around the local minimum With the assumption that the vibrational and the three rotational modes of the 
whole system are decoupled, the full energy of the state s is, to second order in the displacements, £^s,i = f i — rf, 

i—1 K. ) i—Z ' 
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Fig. 8 (color online) Experimental and theoretical probability of stationary states for a Yukawa ball with A'^ = 31. Left: 
Experimental results including mean probability (full lines) and standard deviations (dashes) corresponding to the states (from 
top to bottom) with 5, 4 and 6 particles on the inner shell compared to LMD results (lines with symbols) versus friction at 
K = 0.8. Right: Analytical theory (lines) compared to MC results (symbols) as a function of temperature. The solid lines are 
obtained from Eq. ijTSj and the dashed lines neglect the statistical weight factor caused by the eigenfrequencies, (Qs = 1) for 
all states. Agreement with the experimental data is obtained for T ~ 0.025. 



The term in parentheses denotes the vibrational kinetic and potential energy of the / = 3A^ — 3 normal modes 
with the normal coordinates Even though the harmonic approximation is only applicable for low temperatures 
(or strong coupling) when the particles oscillate around the equilibrium positions with a small amplitude, it can 
be used to estimate phase transitions and compare the results with simulations and experimental data. With the 
factorized partition function 

Zs = nsZl'"*Z^'''Zl°\ with the degeneracy factor n, = , (13) 

Oi— 1 ■ 

where L is the number of shells and the occupation number of shell i with X^^Li — ^- The degeneracy 
factor Us denotes the number of possibilities to form a configuration with shell occupation (iVi, N2, N^) from 
N distinguishable particles. The internal partition function is given by Z™* — e~^^^ , and the vibrational and 
rotational partition functions are 

where we introduced the geometric mean eigenfrequency fls = {Y[i=i '^s,i)^^'^ and the mean moment of inertia 
Is = {Is.ils.2ls.3y^^- Then the occurence probability of the stationary state s is obtained as 

„ 7Xq Z t; Ti e Z _ 

Ps = ^= \ ■ (15) 

The moments of inertia cancel to a good approximation for clusters with 27-40 particles, while for lower particle 
numbers they should be included. Since the total partition function Z is not known, it is advantageous to compute 



instead the ratio of probabilites of two states s and s' , using Eqs. ( 14 1 and ( 15 1 



Thus the probability ratio of two states depends on three factors: their energy difference — , the ratio of 
degeneracy factors Ug/iis' and the ratio of mean eigenfrequencies il^' /^s- The Boltzmann factor e"'^*^^""^^'-' 
gives preference to states with a low energy. For low temperatures it will be the dominant factor but it becomes 
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less important for higher temperatures when fc^T >> E'^, — and e^^^^"^-^"'") « 1. The degeneracy factor 
assigns large statistical weight to states with more particles on inner shells. For example, for = 31 the ratio 
of the factors for two configurations with, respectively, 4 and 5 particles on the inner shell is n(27,4)/n(26,5) = 
(26!5!)/(27!4!) = 5/27. The energy differences are computed by MD simulations and, as an exmmple for 

= 31 and k — 1.4, are given in table [T] The effect of the mean eigenfrequency, i.e. the effect of the local 
curvature of the potential surface cannot be that easily predicted but the ratio of the eigenfrequencies are in the 
range of 0.8 — 30 and can therefore strongly influence the ratio of the probabilty of states. The analytical results 
can be compared with MC simulation results [90] which are practically exact. In particular, they include all 
anharmonic corrections. For a fixed temperature the occurrence probability of each stationary state is recorded. 
Both results are shown in Fig.[8]as a function of temperature, for the example N — 31. Evidently, the analytical 
results are accurate at temperatures below fc^T < 0.02. Above this value anharmonic corrections become 
important. For other particle numbers and screening parameters simliar results are obtained. 

Let us now compare to the experimental results, cf. left part of Fig. |8] These data were obtained at room 
temperature, corresponding to approximately T ~ 0.0015. At this temperature the MC and analytical results 
both are in agreement with each other (right figure) predict practically 100 percent occurence of the ground 
state (27, 4) which is in striking contrast to the experiment. To resolve this contradiction, independent Langevin 
dynamics simulations have been performed [90, 113] which model the experimental situation of turning on and 
off the confinement by including a time-dependent potential [7°''* in Eq. ([T]l and which are included in the left 
figure part. Using a screening parameter which corresponds to the experimental parameters (the same as in the 
MC simulations) one finds a strong dependence on the chosen values of the friction coefficient ly. For small 
the simulation always favor the ground state which is easily understood: with small friction the system is cooled 
very slowly and has enough time to find the lowest potential minimum (this is nothing but simulated annealing 
which was applied in Sec. 4.1 to obtain the ground states). In contrast, with large damping the particles rapidly 
loose their kinetic energy even before settling in a potential minimum. The experimental damping parameters 
are well above i/ = 2 where Fig.[8]indicates very good agreement between LMD and experiment. On the other 
hand, Monte Carlo and analytical theory are based on a hamiltonian model, cf. Sec. |3.2| and do not include 
damping effects. Thus, a direct comparison with the experimental probabilities is not possible. However, as was 
demonstrated in Ref. [90] due to the large friction in the experiments the plasma particles rapidly relax towards a 
Maxwellian distribution already at temperatures of about T = 0.05. During further cooling the system remains 
in equilibrium, only its temperature is being gradually reduced until it reaches T ^ 0.002. The "decision" about 
what stationary state will be reached is "made" at an early stage, before the Yukawa ball freezes which occurs 
around T ~ 0.02. In fact, Fig.[8]shows that the MC simulations at T « 0.025 agree very well with the experiment 
and the LMD [90]. 



6.3 Structural transitions in the experiment 

Transitions between two stationary configurations which were discussed above could, in fact, be observed ex- 
perimentally [113]. Figure|9|a) shows the motion of a single particle of a Yukawa ball with 31 particles over a 
period of ten minutes. It can be seen that the particle moves from the outer shell inwards and returns to the outer 
shell within about 200 seconds. This transition is temperature induced (it occurs spontaneously) and is observed 
in other Yukawa balls as well [113]. This behavior is well explained by the MC simulation of the total potential 
profile [see the discussion in Sec. |6.1[ which are shown in the right part of Fig. |7] The potential barrier for 
moving a particle inwards from (27, 4) to (26, 5) [lowest curve] shows where the particle returns; interestingly 
the sharp increase at tq w 0.85 is matching the return point of r = 0.39mTO in the experiment, cf. Fig.|9] This 
peak corresponds to a MC run with fewer steps where the other particles have not enough time to fuUy relax in 
response to the moving paricle - this nonequilibrium situation is apparantly observed in the experiment. 



6.4 Finite temperature effects. Solid-liquid phase transition 

So far we have considered only the low-temperature behavior of the Yukawa balls where the system was in the 
ground state or low-lying metastable state. Finite thermal energy caused only occasional transitions between 
stationary states. This behavior is typical for matter at very strong coupling characterized by very low kinetic 
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Fig. 9 (a) Radial component of the trajectory of one particle leaving the outer shell of a Yukawa ball with N = 31, the 
inset shows the trajectory in the p-z plane (black symbols), (b) Average Wigner-Seitz cell volume per particle, Vws, for 
all particles belonging to the inner shell as a function of the radial position of the traveling particle. For R = 0.55mm, the 
Wigner-Seitz cell of the traveling particle is included in the average. The solid lines are the average Vws of the configurations 
(4, 27) and (5, 26) obtained by repeating the experiment and averaging. 



energy, typical for a Wigner crystal. This physical situation can be quantified by introducing a coupling param- 
eter, which is defined as the ratio of the mean interaction energy to kinetic energy. For a classical macroscopic 
Coulomb system in thermal equilibrium the definition is 

r - 17 

^ AneaksT' 

where a is the mean interparticle distance. The transition from strong coupling where the system is in a crystalline 
state to a liquid occurs at critical value Tcr which is of the order 175 (137) in 3D {2D). The transition to a gas-like 
state occurs around Vc — 1' (for generalizations to quantum systems, see Ref. [117]). 

Here we are concerned with systems with Yukawa interaction and have to generalize the coupling parameter 
The straighforwad way to do this would be to replace the Coulomb energy of a pair of particles by the mean 
Yukawa energy and Tc e^*^" • Fp, but this makes it impossible to systematically define a value for the melting 
transition in a Yukawa system, so for the phase diagram [118] one has to use other quantities. On the other hand, 
it is not the total potential energy which governs the melting transition, but the process is driven by the excitation 
energy of particles above the potential minimum which leads to increasing fluctuations of the particles around 
their stationary positions. These fluctuations are defined by the local curvature of the potential [99, 117], i.e. its 
second derivative in the minimum. This allows to define a proper effective Yukawa coupling parameter [1 19], 

which allows to qualitatively correctly predict the melting point Vcr{i^) for any k, based on the known critical 



value for a Coulomb system, i.e. by using, in Eq. ( 18 i, Tc — Tcr- 

However, for the analysis of melting in Yukawa balls one has to take into account the peculiarity of small 
systems. Finite systems behave fundamentally different to macrososcopic systems: there exist no strict phase 
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transitions. Nevertheless, a similar transition or crossover from solid-like to liquid-like behavior is well known 
and has been studied as a function of particle number e.g. by SCHIFFER [95]. For small N the melting tem- 
peratures are typically higher than in the macroscopic limit and the crossover proceeds over a finite range of 
coupling parameters (temperatures) making it often difficult to identify a melting point, see below. Furthermore, 
finite systems are strongly influenced by the symmetry of the confinement potential, and the shell structure of the 
Coulomb balls affects the melting behavior 

The effect of screening on the stability and melting of Yukawa balls was analyzed in Ref. [99], where a simple 
analytical estimate can be observed for the example N ~ 2. Expanding the total potential energy to second order 



around the ground state, as in the macroscopic system, a formula analogous to ( 18 i is obtained 



"(-'^-^)-"-(^-^) l + ..o. + .%V/3 - '''' 
Here rgy (k) is the two-particle distance for Yukawa interaction which is obtained by inverting the equation 

l:r^ = i^-rUroy), (20) 
1 + KroY 47rem 

where vq is the equilibrium distance of two particles with Coulomb interaction, see Sec. [3] As in the macroscopic 
case, screening always destabilizes the crystal. 

Another quantity to locate the melting point are the mean position fluctuations of particles which are known 
to be small (large) in the solid (liquid) phase, independently of the interaction and system size [117] which is 
commonly known as Lindemann criterion. However, this parameter diverges in 2D systems, and an improved 
quantity are the mean relative interparticle distance fluctuations, i.e. the fluctuations normalized to the mean 
distance 




This quantity has been successully applied to many finite systems and also to Yukawa balls [89,98, 100, 107], an 



example is shown in the left part of Fig. 10 One clearly sees the increase of the fluctuations typical for a melting 
process. It is also evident that this increase extends over more than one order of magnitude in the temperature. 
A rather detailed analysis of melting in Coulomb balls was performed by APOLINARIO [100] who observed that 
the melting often proceeds in several steps. 

Finally, in a recent study [120] it was observed that the quantity Urei has a number of difficulties in small 
systems - the result strongly depends on the way of calculation. Instead it was suggested to use another quantity 
- the variance of the block averaged interparticle distance fluctuations (VIDF) 

1 ^ I 

^u..^ = ^ E \l^<eli^)) - i^relis))^. (22) 
s=l 

This parameter measures the width of the probability distribution P{urei) which has a maximum in the vicinity 
of the phase transition. In the solid phase the particles show only small distance fluctuations, and in the liquid 
phase all particles strongly fluctuate and P{urei) has a peak centered around a large value, yet the width of the 
peak is again small. However, in the vicinity of the phase transition, the particles may be in one state for some 
time and frequently make transition over the potential barrier to another stable state and so on. This causes a 
broad peak of P, as can be seen in the rightmost part of Fig.[TO| 

An experimental analysis of the melting behavior is difficult because it require a systematic scan of the plasma 
parameters in a sufficiently broad range around the expected melting point. To do this without changing the 
screening is difficult and requires future advances in the handling of the Yukawa balls. 

6.5 Further results on Yukawa balls 

In this review we had to concentrate on a few properties of Yukawa balls leaving out many other interesting 
results some of which are mentioned in the following. 
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Fig. 10 Mean distance fluctuations IzTl, crosses and their variance i22 ; circles, for a Yukawa ball with k — 0.4 in the 
vicinity of the melting point Tcr- Left: i'emperature dependence of Urei and (7„^,^j for a Yukawa ball with A'' — 40. Right: 
Block averaged interparticle distance fluctuations Urei (s) vs block number s during an MC simulation and their accumulated 
probability P (rightmost figure) for N = 40 and four temperatures (from bottom to top): below, close to, at and above Tcr- 
Block length M = 1000. 



The normal mode spectram of Yukawa balls which was mentioned in Sec. 6.2 has been analyzed in great 
detail. Among the most interesting collective excitations of these finite systems is the breathing mode (uniform 
radial expansion and contraction). Recently it could be shown [121] that Yukawa balls do not possess such a 
mode. Instead they can exhibit simultaneously several monopole oscillations (periodic increase and reduction of 
the mean radius). Also, first experimental observations of the normal modes of small Yukawa balls have been 
performed [122]. 

The confinement potential in the experiments is certainly never perfectly isotropic, so the analysis of anisotropy 
effects is of interest. Ground and metastable states in anisotropic parabolic confinements have been analyzed as 
a function of the anisotropy by APOLINARIO et al. [98, 100]. They also analyzed the normal modes of these 
systems. The results are of interest also for experiments where the confinement is made anisotropic on purpose. 

Another interesting theoretical investigation has been performed by PSAKHIE et al. They performed MD 
simulations for Yukawa balls containing two sorts of particles having different masses [97]. In their model 
the confinement forces, the termophoretic force, gravity, electric field force and ion drag have been explicitely 
included. Increasing the mass difference they found that beyond a critical mass asymmetry the Yukawa ball may 
split in two balls containing different particles which are vertically separated. 

Recently, first simulations of the short-time dynamics of Yukawa balls have been performed by Kahlert et 
al. [90]. Starting from a random initial state and rapidly turning on the confinement showed a very interesting 
nonequilibrium behavior with a non-monotonic kinetic energy increase. This behavior is very different from the 
ultrafast dynamics of weakly coupled plasmas and is a very interesting subject of further research. 



7 Discussion and Outlook 

In this review, and overview on the recently discovered spherical crystals in dusty plasmas - Yukawa balls - has 
been given. The experimental setup and the diagnostic methods have been briefly discussed. The main attention 
was then devoted to the theoretical modeling and to computer simulatinos of Yukawa balls, when possible, aiming 
at a detailed comparison to the experimental data. This close connection between experiment and theory has lead 
to substantial advances in our understanding of finite Yukawa crystals. It was shown that the Yukawa balls 
produced in the Kiel and Greifswald labs are well described by an isotropic statically screened pair potential and 
that wake effects are apparantly not of importance. This situation may change when larger balls are studied which 
is an interesting question for further research. 
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It was shown that screening has a crucial effect on the internal structure of Yukawa balls, including the shell 
populations, the ground and metastable states. Comparisons between measurement and simulations allowed to 
accurately determine the value of the screening parameter. Besides the shell structure, also the mean distribution 
of charge throughout the Yukawa balls could be explained. This was achieved by an analytical theory based 
on a fluid-like description. It could be shown that, with increased screening, the average density is more and 
more concentrated in the central region of the cluster and rapidly decays towards the edge. This is an important 
difference to Coulomb systems, including spherical ion crystals, which exhibit a constant radial density profile. 

Not only the grovmd state of the Yukawa balls could be observed. The experiments were able to reUably 
detect metastable states and to reveal the probability of their occurence. These measurements could be well 
explained by theory and simulations where it became clear that the (in some cases) dominant role of metastable 
states - despite relatively low temperature - is due to the large friction experienced by the dust particles in the 
plasma. Besides the static properties also several dynamic aspects are now well understood. The experiments 
could demonstrate that spontaneous transitions between various stationary states occur which is well explained 
by a theoretical analysis of the potential energy landscape, in particular from the energy barriers between local 
potential minima. The theory further predicted a rather interesting melting behavior of the Yukawa balls upon 
further temperature increase. It will be of high interest to verify these predictions in future experiments. 

At the same time many fundamental theoretical questions remain to be answered. In the shell strucutre of 
the Yukawa balls correlation effects play a crucial role. This became clear by studying greatly simplified shell 
models. These effect were also important to explain the density profile of Yukawa balls at large screening using 
a local density approximation. Yet no analytical results for the pair correlations are known and no derivation 
of these effects from statistical mechanics has been given. The excellent experimental accuracy gives a unique 
opportonity to advance and verify modern theoretical approaches in strongly coupled plasmas. 

FinaUy, questions of further interest in the field of Yukawa balls include their normal mode spectrum, and their 
interaction with external electric and magnetic fields. This will require to take into account the modification of 
the pair interaction. Here also wake effects could play a more important role. Last but not least it is becoming 
possible to advance to an analysis of the short-time behavior Yukawa balls. Excerting a fast excitation, e.g. 
by laser manipulation, the time-dependent response of these strongly correlated systems can be systematicaUy 
studied. 
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